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Abstract 


MERRA products were used to force an established ocean biogeochemical model to estimate 
surface carbon inventories and fluxes in the global oceans. The results were compared to public 
archives of in situ carbon data and estimates. The model exhibited skill for ocean dissolved 
inorganic carbon (DIC), partial pressure of ocean CO 2 (PCO 2 ) and air-sea fluxes (FCO 2 ). The 
MERRA-forced model produced global mean differences of 0.02% (approximately 0.3 pM) for 
DIC, -0.3% (about -1.2 patm; model lower) for pCC> 2 , and -2.3% (-0.003 mol C m~ 2 y' 1 ) for 
FCO 2 compared to in situ estimates. Basin-scale distributions were significantly correlated with 
observations for all three variables (r=0.97, 0.76, and 0.73, P<0.05, respectively for DIC, pCC> 2 , 
and FCO 2 ). All major oceanographic basins were represented as sources to the atmosphere or 
sinks in agreement with in situ estimates. However, there were substantial basin-scale and local 
departures. 
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1 Introduction 


The Modem-Era Retrospective analysis for Research and Applications (MERRA) project 
represents a next generation of reanalysis products. Utilizing data from NASA Earth observing 
satellites, MERRA is intended to improve upon the widely recognized set of existing reanalysis 
products, primarily by including a more realistic representation of the hydrological cycle 
(Rienecker et ah, 2011). A comprehensive approach using advanced data assimilation 
methodologies and modem Earth remote sensing observations, along with state of the art 
atmospheric and hydrological models, MERRA is expected to fully support climate-related 
modeling efforts. 

Here we use MERRA products to force a global ocean biogeochemical model. We seek to 1) 
simulate the distributions and fluxes of carbon components in the global oceans with an explicit, 
prognostic description of the carbon cycle and 2) evaluate the realism of the model results. The 
simulation is accomplished using an established three-dimensional model of the global oceans 
containing prognostic representations of biological and chemical constituents involved in the 
ocean carbon cycle. Evaluation is achieved through comparison with observations of surface 
carbon inventories and fluxes. 


5 


2 Methods 


2.1 Global Three-Dimensional Circulation Model 

Global ocean carbon dynamics are simulated by the NASA Ocean Biogeochemical Model 
(NOBM; Figure 1). A complete description of the model can be found in Appendix 1. It is a 
three-dimensional representation of coupled circulation/ biogeochemical/radiative processes in 
the global oceans (Gregg et al., 2003; Gregg and Casey, 2007). It spans the domain from -84° to 
72° latitude in increments of 1.25° longitude by 2/3° latitude, including only open ocean areas, 
where bottom depth>200m. The biogeochemical processes model contains 4 phytoplankton 
groups, 4 nutrient groups, a single herbivore group, and 3 detrital pools (Figure 2). The 
phytoplankton groups differ in maximum growth rates, sinking rates, nutrient requirements, and 
optical properties. The 4 nutrients are nitrate, regenerated ammonium, silica to regulate diatom 
growth, and iron. Three detrital pools provide for storage of organic material, sinking, and 
eventual remineralization. 


NASA Ocean Biogeochemical Model (NOBM) 
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Figure 1. Interactions among the main components of NOBM, nominal outputs, and forcing fields. IOP indicates 
inherent optical properties. Forcing variables are shown in the gray boxes. MERRA forcing variables are in bold. 
Surface pressure and precipitable water effects on surface irradiance play a small role in the inorganic carbon results 

and are ignored in this effort. 

Carbon cycling involves dissolved organic carbon (DOC) and dissolved inorganic carbon 
(DIC; Figure 2). DOC has sources from phytoplankton, herbivores, and carbon detritus, and a 
sink to DIC. DIC has sources from phytoplankton, herbivores, carbon detritus, and DOC, and 
communicates with the atmosphere, which can be either a source or sink. The ecosystem sink 
for DIC is phytoplankton, through photosynthesis. This represents the biological pump portion 
of the carbon dynamics. The solubility pump portion is represented by the interactions among 
temperature, alkalinity (parameterized as a function of salinity), silica, and phosphate 
(parameterized as a function of nitrate). The alkalinity/salinity parameterization utilizes the 
spatial variability of salinity in the model adjusted to mean alkalinity 
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TA = TA S/S 


where TA is total alkalinity and S is salinity. The underscore represents global mean values. TA 
is specified as 2310 pE kg' 1 (Ocean Model Intercomparison Project (OCMIP; 
www.ipsl.jussieu.fr/OCMIP) and S as 34.8 PSU (global model mean). Since the model contains 
nitrate but not phosphate, we adjust nitrate by multiplying by 0.1. This is derived from the 
global mean ratio of nitrate to phosphate from the National Oceanographic Data Center 
(Conkright et al., 2002) for their top three standard levels. The calculations for the solubility 
pump follow the standards set by the Ocean Model Intercomparison Project (reference above). 
We employ a locally-developed lookup table valid over modem ranges of DIC, salinity, 
temperature, and nutrients for computational efficiency, at no cost to accuracy. Air-sea CO 2 
exchange uses the Wanninkhof (1992) formulation, as is common in global and regional ocean 
carbon models (e.g., McKinley et al., 2006). 



Figure 2. Pathways and interactions among the components of the biogeochemical processes model, illustrating the 
interactions with the carbon cycle, comprising dissolved inorganic carbon (DIC), dissolved organic carbon and 
exchanges with the atmosphere as a function of the ocean and atmosphere partial pressures of C02 (pC02). The 
biological pump is represented by phytoplankton, herbivores, nutrients, and detritus. 

We acknowledge that the parameterization of alkalinity using salinity is a simplification. 
However, the relationship in the modem oceans is robust, especially in tropical and subtropical 
oceans (Millero et al., 1998; Key et al., 2004; Lee et al., 2006) and has been recommended for 
international carbon model intercomparison efforts (www.ipsl.jussieu.fr/OCMIP). Furthermore, 
the salinity-parameterization of alkalinity used in the model compares favorably globally and 
over major oceanographic basins with in situ data sets for the modem oceans (Figure 3). The 
basins are shown in Figure 4. 

NOBM undergoes spin-up for 200 years under climatological forcing from MERRA. Initial 
conditions for DIC are derived from the Global Data Analysis Project (GLODAP; Key et al., 
2004). We average DIC over oceanographic basins and depth and use these mean values for 
initial conditions. DOC initial conditions are set to 0 pM. Other initial conditions are described 
in Gregg and Casey (2007). The first ten years of the mn show a net pC02 difference ApC02 
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(year 10-year 1) of -0.982 juatin, at the first hundred years the 10-year ApCC >2 (year 100-year 91) 
is 0.413 patm, and at 200 years, the 10-year ApCC >2 (year 200-year 191) is 0.102 patm. 
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Figure 3. Annual surface alkalinity distributions and statistics from the model and data. The model-data basin 
correlation is statistically significant at P<0.05, as indicated by an asterisk. 



Figure 4. Delineation of the 12 major oceanographic basins. 
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2.2 Data Sets 


2.2.1 Forcing Data Sets 

Forcing data sets are shown in Figure 1. Monthly climatologies are used in all cases. All 
except soil dust (iron), ozone, clouds, and atmospheric CO 2 are obtained from MERRA 
products. Ozone is from the Total Ozone Mapping Spectrometer and Ozone Monitoring 
Instrument, and soil dust deposition is from Ginoux et al. (2001). Cloud data (cover and liquid 
water path) are obtained from the International Satellite Cloud Climatology Project. 
Atmospheric CO 2 is taken from the Lamont-Doherty Earth Observatory (LDEO) data set 
(Takahashi et al., 2009), using a mean over the entire range of observations of 358.7 patm. 
Although the ocean pC02 observations are nominally normalized to the year 2000 (Takahashi et 
al., 2009), we keep the mean atmospheric value from the data to help represent the spatial 
variability at the time of measurement. Our emphasis is to understand how the simulated surface 
carbon components compare with in situ data sets, so using atmospheric values from these 
datasets facilitates the comparison. 

Comparison Data Sets 

The main outputs of interest in this effort are dissolved inorganic carbon (DIC), partial pressure 
of CO 2 (PCO 2 ), and the flux of CO 2 (FCO 2 , notation following Arrigo et al., 2010; Doney et al., 
2009), representing the exchange of carbon between the atmosphere and ocean (positive is 
defined upward, indicating a source to the atmosphere). DIC data sets are obtained from 
GLODAP (Key et al., 2004), which are mapped on a 1° horizontal grid with 33 standard depth 
levels (http://cdiac.oml.gov/oceans/glodap/). pCC >2 and FCO 2 data sets are mapped on a 5° 
longitude by 4° latitude horizontal grid and are surface only. They are obtained from the 
Lamont-Doherty Earth Observatory (LDEO) 

(http://cdiac.oml.gov/oceans/LDEO_Underway_Database/index.html; Takahashi et al., 2009). 


2.2.2 Methodological Approach 

The global model NOBM is forced with the MERRA variables shown in bold in Figure 1. The 
model is spun up for 199 years using monthly climatologies of MERRA forcing. In the 200 th 
year of simulation, model results of surface ocean carbon are compared graphically and 
statistically with climatological in situ data sets and estimates from GLODAP (DIC) and LDEO 
(pC0 2 and FCO 2 ). Results are evaluated globally and regionally in 12 major oceanographic 
basins (Figure 4). Statistical comparisons include global and basin differences between model 
and data global and regional means, expressed as percent, and correlation analysis. Our 
emphasis is on large scale results, so our correlation analysis is performed across the basins 
(N=12, with 10 degrees of freedom). All analyses here are performed for annual mean results, 
and the data sets are converted to the NOBM spatial grid prior to comparison. 
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3 Results 


Surface DIC from NOBM compares favorably with in situ data (Figure 5). There is substantial 
geographical similarity, and basins follow similar patterns. The global difference is 0.02%, 
representing ADIC=0.3|uM, with a correlation across basins of 0.97, which is statistically 
significant at P<0.05. The most notable differences are the tropical upwelling regions. 
Upwelling in the Equatorial Pacific is not apparent in the in situ data, while it is prominent in the 
model. The in situ data show a depression in DIC in the Equatorial Atlantic upwelling region, 
which is in contrast to an increase seen in the model. The basin mean, however, is in agreement 
(Figure 5). The model also has modestly lower DIC in most of the Antarctic. The largest basin 
difference is 1.2% in the Equatorial Indian and Pacific (ADIC=24.5pM and ADIC=23.3pM, 
respectively). The remaining basin differences are <1%. 
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Figure 5. Annual surface dissolved inorganic carbon (DIC) distributions and statistics from the model (blue) and 
data (green). The model-data basin correlation is statistically significant at P<0.05, as indicated by an asterisk. 
Southern, tropical, and northern basins are delineated on the bar graph. Error bars indicate the standard deviation. 
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There is also considerable similarity between the model and in situ data in global pC02 
distributions (Figure 6). The global difference is -0.3% (model lower), which is a difference of - 
1.2 pat in. The correlation across basins is statistically significant with r = 0.76. The model 
pC0 2 distributions are in general agreement in the tropical Pacific upwelling, in contrast with 
DIC, and there is no depression in the tropical Atlantic in the pC0 2 data. There is an east-west 
departure in the North and Equatorial Indian basins between the model and data, with lower 
values in the model in the west, and higher values in the east. There is much more spatial 
variability in the model, as seen in the maps and also the standard deviations (Figure 6). The 
largest basin difference is -8.9% in the South Atlantic ( ApC O 2 =- 3 2 . 3 p atm) followed by the 
Equatorial Indian at 4.3% (ApC0 2 =16.1patm). The remaining basin differences are <4%. The 
largest local discrepancies are in the eastern South Pacific and Atlantic, where the model exhibits 
two regions of depressed pC0 2 relative to the data, representing about a 35% difference. 
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Figure 6. Annual partial pressure C0 2 (pC02 patm) distributions and statistics from the model (blue) and data 
(green). The model-data basin correlation is statistically significant at P<0.05, as indicated by an asterisk. Southern, 
tropical, and northern basins are delineated on the bar graph. Error bars indicate the standard deviation. 

There are many areas of agreement in the flux of CO 2 (FCO 2 ) between the model and in situ 
estimates, but there are also many areas of disagreement and they are more pronounced than with 
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DIC and pCC >2 (Figure 7). The global difference is -2.3%, representing AFC02=0.003 mol C m" 
2 y' 1 . Like DIC and pCC> 2 , the basin correlation is statistically significant (r=0.73, P<0.05). The 
sign of the flux (source or sink) is the same for model and in situ estimates in all the basins, but 
there are large differences in some basins, such as the South Atlantic and Pacific, North Indian, 
Equatorial Atlantic, and North Pacific. 
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Figure 7. Annual C0 2 flux (FC0 2 ) distributions and statistics from the model (blue) and data (green). The model- 
data basin correlation is statistically significant at P<0.05, as indicated by the asterisk. Southern, tropical, and 
northern basins are delineated on the bar graph. A positive flux indicates a source to the atmosphere. Error bars 
indicate the standard deviation. 


There is much greater spatial variability in the model than the in situ estimates of air-sea flux, 
as indicated in the maps and the standard deviations of the basin means (Figure 7). Locally, 
there are major sinks represented in the model in the northeastern portions of the South Pacific 
and Atlantic that do not appear in the data estimates (Figure 7). The Antarctic exhibits regions of 
strong sources that are also not apparent in the data. The location of sources and sinks in the 
North and Equatorial Indian basins are switched in the model from the data, where the model 
shows a sink in the west and a source in the east, while the data show the reverse pattern. 
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4 Discussion 


The MERRA- forced biogeochemical model (NOBM) produces global estimates of DIC, pCC> 2 , 
and FCO 2 to within a few percent of corresponding in situ data sets and estimates. The global 
mean difference is 0.02%, -0.3%, and -2.3% respectively for DIC, PCO 2 , and FCO 2 , and all show 
statistically positive correlation with data (P<0.05) across the 12 major oceanographic basins. 
Our emphasis here is global and basin scales representations of ocean carbon surface inventories 
and fluxes, as appropriate for a global model, and the results suggest that NOBM, forced by 
MERRA atmospheric and oceanic variables, possesses skill for simulating ocean carbon 
quantities at these large spatial scales. 

We acknowledge that the agreement for DIC is perhaps not surprising, since it was used for 
initial conditions. But it is noteworthy that the model is able to hold these large scale 
distributions over 200 years, and also exhibit substantial agreement with smaller scale 
distributions (Figure 5). We note that, in contrast to pCC >2 and FCO 2 , the spatial variability of 
model DIC, as represented by the standard deviation, is similar to in situ data on global and basin 
scales. 

The global scale agreement of carbon estimates from the model with in situ sources is 
encouraging, and is the main emphasis here, but it is notable that at smaller scales the agreement 
is not as good. Basin and local discrepancies increase as one evaluates the carbon from DIC to 
pCC >2 to FCO 2 . The increase in basin and local scales from pCC >2 to FCO 2 is a reflection of the 
fact that small discrepancies between model and data in DIC and pC02 can produce important 
differences in FCO 2 . Thus FCO 2 is sensitive to small errors in DIC and pCC> 2 , and it can be a 
challenge to represent this well in models, especially at smaller scales and using a global model. 

We note that use of atmospheric pCC >2 for year 2000 only modestly changes the results shown 
here. Specifically, the global mean difference is 0.3%, 2.2%, and 7.8% respectively for DIC, 
pCC> 2 , and FCO 2 , and again all show statistically positive correlation with data (P<0.05) across 
the 12 major oceanographic basins. The differences represent a stronger flux into the oceans as 
expected from the higher atmospheric pCC >2 (368.6 instead of 358.7 patm). 
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6 Appendix A 

6.1 NASA Ocean Biogeochemical Model - Equations 

The governing equations of the model are 

Phytoplankton 

3 3 

-P t = V-CKV^-WPi-— (wJjPjT Qq — (6 + O)] Pj — yH — kP s (Al) 


i = 1 => diatoms 
i = 2 => chlorophytes 
i = 3 => cyanobacteria 
i = 4 => coccolithophores 

Nutrients 

-N h = VCKVN^-VVN^-bf,^ + Ra c D c / (C: N) + A d D c / (C: N) 

3 

- N A = V-(KVN a ) - VVN a - b N [X, PiPj + b N s[yH + n 2 H 2 ] 

^N s = V.(KVN s )-V.VN s -b s ji 1 P 1 + Ra 3 D 3 
3 

-N F = V(KVN F )-V*VN F -b F [X i p i P i ] + bp £ [yH + fi,H 2 ]+Ra F D F + A Fe /L-9N F 

at 

N n = nitrate 
Na = ammonium 
N s = silica 
N F = dissolved iron 

Herbivores 

3 

— H = V»(KVH) - VVH - (1 — s)yH — r^H — r| 2 H z - toH - 0H 


(A2) 

(A3) 

(A4) 

(A5) 


(A6) 


Detritus 

3 9 

-D c = V.(KVD c )-V.VD c - — (w d ) c D c -Ra N D c + * [kE* P t + i^H] + *(1 - e^H 2 


dt 


3z 


(A7) 


^-D^C 
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d d 

-D s = V(KVD S > - VVD S - — (w d ) s D s - Rot s D s + b s [kP 1 + yH] 

3 d 

- D F = V»(KVD f ) - V«VD F - — (w d ) F D F - Ra F D F + b F [kE l P L + 1 ] iH] + b F (1 - €)r| , H 2 
ot o z 

+ 0Np 

Dc = carbon/nitrogen detritus 
D s = silica detritus 
D F = iron detritus 


(A8) 

(A9) 


Carbon 


^ DOC = V»(KVDOC) - ^VDOC + O d^P; + <D G)H + k D D c - cpDOC (A10) 

d d Ra N D c 

— DIG = V(KVDIC)- — VDIC-4>2H i P i +^EM i P i + ^0H + «pDOC+— ^ + AO co ^ (All) 
£7t az (C:INJ 


where the symbols and values are identified in Table 1. Bold denotes a vector quantity. All 
biological processes are assumed to cease in the presence of sea ice, in proportion to the fraction 
of sea ice cover, which is included as an external forcing field. 


6.2 NASA Ocean Biogeochemical Model - Description 

Circulation Model 

The Ocean General Circulation Model (OGCM) is a reduced gravity representation of circulation 
fields (Schopf and Loughe, 1995). It is global in scale, extending from near the South Pole to 
72° N, in increments of 2/3° latitude and 1 1/4° longitude, comprising all regions where bottom 
depth > 200m. The model contains 14 vertical layers, in quasi-isopycnal coordinates, and is 
driven by wind stress, sea surface temperature, and shortwave radiation. 

Radiative Model 

Radiative transfer calculations provide the underwater irradiance fields necessary to drive 
growth of the phytoplankton groups, and interact with the heat budget. The Ocean-Atmosphere 
Radiative Model (OARM; Gregg, 2002) contains a treatment of the spectral and directional 
properties of radiative transfer in the oceans, and explicitly accounts for clouds. The 
atmospheric radiative model is based on the Gregg and Carder (1990) spectral model, extended 
to the spectral regions 200 nm to 4 pm. It requires external monthly climatologies of cloud 
properties (cloud cover and liquid water path), surface pressure, wind speeds, relative humidity, 
precipitable water, and ozone. Aerosols are considered to be strictly of marine origin and are 
computed as in Gregg and Carder (1990). 
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Oceanic radiative properties are driven by water absorption and scattering, the optical 
properties of the phytoplankton groups, and chromophoric dissolved organic matter (CDOM). 
Three irradiance paths are enabled: a downwelling direct path, a downwelling diffuse (scattered) 
path, and an upwelling diffuse path. All oceanic radiative calculations include the spectral nature 
of the irradiance. 

Optical properties of coccolithophores and other phytoplankton groups were derived from 
laboratory studies. Their values and references can be found in Gregg (2002). 

Phytoplankton 

The growth formulation includes dependence on total irradiance (E T ), nitrogen as nitrate plus 
ammonium (N T ), silica (Si - for diatoms only), iron (Fe), and temperature (T) 

Pi = (fW)i min [f(N T ,Si,Fe, £V) ( ] RG i , (A12) 


where i indicates the phytoplankton functional group index (in order, diatoms, chlorophytes, 
cyanobacteria, and coccolithophores), p is the total specific growth rate (d' 1 ) of phytoplankton, 
p max is the maximum growth rate at 20°C (Table 1). The term p(E T ) represents the growth rate 
as a function solely of the total irradiance (pmol quanta m' s' ), 


fi(E T ) = 


l f 


(A 13) 


where k E is the irradiance at which p = 0.5p m and equals 0.5 h, where E is the light saturation 
parameter. The nutrient-dependent growth terms are 


fi(N0 3 ) = 


NO, 


(A 14) 


fiiNVj = 


(NH 4 + (k N ) t ) 


(A 15) 


f t (N r ) = fi(NH 4 ) + min [/;(NO g ),l] - 


(A 16) 


(Gregg and Walsh, 1992) 




Si 


[5i + (fcj J 


(A 17) 


fi( Fe ) 


Fb 

[Fs-l-Ck F } j ’ 


(A 18) 


Temperature-dependent growth is from Eppley (1972) 

R = 1 . 066 ^ T_2trc ^ 


(A 19) 
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which produces a temperature-growth factor normalized to 20°C. The term G in Eq. A12 is an 
additional adjustment used for the cyanobacteria component that reduces their growth rate in 
cold water (<15°C) 


G, = 0.0294T + 0.558 . 


(A20) 


G; = 1 for the other three phytoplankton components (i=l,2,4) where 0.0294 is in units of per 
degree Celsius. When T>15°C, G 3 reaches its maximum value of 1. This effect conforms to 
observations that cyanobacteria are scarce in cold waters (Agawin et al., 2000; 1998). 

The cyanobacteria component possesses a modest ability to fix nitrogen from the water 
column, as observed in Trichodesmium spp. (Carpenter & Romans 1991). The nitrogen fixation 
is expressed as additional growth occurring when nitrogen availability is <(k N ) 3 , 

= °' 25 (~ 7SP a) , (A21) 


where the index 3 indicates cyanobacteria. The biomass dependence represents a progressive 
community changeover from non-N-fixing cyanobacteria to N-fixing bacteria as the total 
population declines under nitrogen-stressed conditions. The total N-limited growth rate plus the 
additional growth derived from N-fixation is not allowed to exceed the growth rate where total 
nitrogen = (k N ) 3 . No accounting for denitrification is made in the model. 

Photoadaptation is simulated by stipulating 3 states: 50, 150 and 200 (pmol quanta m' s' ). 
This is based on laboratory studies which typically divided experiments into low, medium, and 
high classes of light adaptation. Carbon: chlorophyll ratios (<f>) are kept constant at 50 g:g for 
mass conservation. 

Mean irradiance is computed during daylight hours, and then the phytoplankton photoadaptive 
state is classified accordingly. This calculation is only performed once per day to simulate a 
delayed photoadaptation response. Light saturation constants for the three light levels are 
provided in Table 1 . 

Phytoplankton group physiological parameters p m , Ik, and k N ,s,F are derived from carefully 
controlled, inter-comparative laboratory studies. We require that at least two of the groups are 
involved simultaneously in order to utilize the experimental results. For p m mean values of the 
relative growth rates are derived from the results of Ben-Amotz and Gilboa (1980), Brand et al. 
(1986, 1983), Eppley et al. (1969), Falkowski et al. (1985), Furnas (1991), Gavis et al. (1981), 
Goldman and Glibert (1982), Humphrey (1979), Subba Rao (1981), and Sunda and Huntsman 
(1995). 

Light saturation parameters, Ik, are formulated for the three irradiance categories used to define 
photoadaptation. Mean values are summarized from the reports of Barlow and Alberte (1985), 
Bates and Platt (1984), Langdon (1987), Perry et al. (1981), Sakshaug and Andresen (1986), and 
Wyman and Fay (1986). 
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The coccolithophore half-saturation constant for nitrogen (k N ) was observed by Eppley (1969) 
to be one-half the value of diatoms. Cyanobacteria k n is set slightly lower than 
coccolithophores, assuming small particle size leads to improved nutrient uptake efficiency. 
Chlorophyte k N is set at one-third the departure between diatoms and coccolithophores. The 
diatom k N is arbitrarily set to 1 pM. 

Phytoplankton vector sinking is treated as additional advection in the z-direction. Sinking rates 
are specified at 31°C and derived from Stokes Law using representative phytoplankton sizes 
from Ahn et al. (1992), Bricaud and Morel (1986), Bricaud et al. (1983; 1988), Dubinsky and 
Berman (1986), Kirk (1980); Mitchell and Kiefer (1988), Morel (1987), Morel and Bricaud 
(1981), and Sathyendranath et al. (1987), for the individual groups. The rates are adjusted by 
viscosity according to Stokes Law (Csanady, 1986), which is parameterized here by temperature 

w 3 (r) = w 3 (31 D C)[a451 + 0.0178r] . (A22) 


Coccolithophore sinking rates are allowed to vary as a function of growth rate from 0.3 to 1.4 m 
d' 1 based on observations by Lritz and Balch (1996). A linear relationship is assumed 

= 0.752 n 4 (higk) + 0.225 , (A23) 

where w s is the sinking rate of coccolithophores (m d 1 ), p(high) is the highest growth rate 
actually achieved for the previous day, and the subscript 4 represents coccolithophores. 


Nutrients 

The diversity in the processes affecting the four nutrient groups requires elucidation in 4 
separate equations, unlike the phytoplankton. All are taken up by phytoplankton growth, with 
silica subject only to diatom uptake (note the subscripts in Eq. A4 denoting diatoms). Lor three 
of the nutrients, nitrate, silica, and dissolved iron, corresponding detrital pools remineralize to 
return nutrients previously uptaken by phytoplankton. There is no detrital pool for ammonium, 
which is excreted as a function of herbivore grazing, and as a function of higher order ingestion 
of herbivores, represented by the term n 2 H 2 in Eqs. A3, A5, A6, A7, and A9. Dissolved iron 
also has an excretion pathway, but nitrate and silica do not. The nutrient to chlorophyll ratios, 
denoted b in Eqs. A2-A5, are derived from Redfield ratios (Table 1) and the carbon: chlorophyll 
(<f>) ratio is constant. 


b v =®/C:N 
b s = Q/C-.S 
b F = e P/CiFe . 


(A24) 

(A25) 

(A26) 


As in Gregg et al. (2003) dust deposition fields are derived from Ginoux et al. (2001). In this 
model, four dust size fractions are transported, corresponding to clay (smallest) and three increasing 
fractions of silt. The iron content is assumed to vary among the clay and silt fractions as follows: 
clay = 3.5% iron, silt = 1.2% iron (Fung et al., 2000). Iron solubility is assumed at 2% for all 
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fractions, which is toward the low end of current estimates (Fung et al., 2000), but is the same as 
used by Moore et al. (2004). The bottom boundary condition is 0.6 nM (Archer and Johnson, 
2000 ). 

Iron scavenging is implemented in this version of NOBM. It is set at 2.74xl0' 5 d' 1 at low iron 
concentrations (<0.6nM; Moore et al., 2002) and 50 times this rate at higher concentrations. A 
smooth transition is enabled as in Moore et al. (2002). 

6 - 2.74 X 10 ~ 5 N F ,for N F < 0.6 nM (All) 

Q - 2.74 X 10 ~ s N f + 1.37 X 10“ 3 (N f - 0.6), for N F > O.SnM . (A28) 


Herbivores 

Grazing uses an Ivlev formulation (McGillicuddy et al., 1995), 

y(r) = [I - exp(-A vp)] . (A29) 

Rh is the maximum grazing rate at 20° C (y m ) adjusted by temperature 

R h = 0.06 exp(O.ir) + 0.70 . (A30) 


The temperature-dependence for grazing is more linear than that for phytoplankton, reflecting 
the larger size of their overall community. The grazing represents the total loss of phytoplankton 
to herbivores, as indicated by the summation symbol, but is applied to the individual 
phytoplankton functional groups proportionately to their relative abundances. This enables 
herbivore grazing to self-adapt the prevailing phytoplankton community. 

The two loss terms in Eq. A6 represent the death of herbivores (r|iH) and higher order 
heterotrophic losses (r| 2 H 2 ). These formulations and parameters (Table 1) were taken from 
McGillicuddy et al. (1995). 


Detritus 

Three detrital components represent the three major nutrient elements, carbon/nitrogen, silica, 
and iron (Eq. A6-A9). The nitrogen detritus is kept as carbon in the model, but since the C:N 
ratio is constant, it is simple to convert when needed. All are subject to advection, diffusion and 
sinking. Detrital sinking, like phytoplankton sinking, is dependent on viscosity parameterized 
here in terms of temperature, using the same formulation. Remineralization, k, is also 
temperature-dependent, and uses the phytoplankton growth-dependence term R in Eq. A 19. 
Silica contained in the diatom component of phytoplankton is assumed to pass through 
herbivores upon grazing directly into the silica detritus pool. No silica remains in the herbivore 
component at any time. 
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Carbon 


Dissolved organic carbon cycling is taken from Aumont et al. (2002), with conversions added 
for compatibility with NOBM units. In addition, all parameters are temperature-dependent, 
unlike Aumont et al. (2002), using the phytoplankton temperature dependence defined in Eq. 

A 19. Following Aumont et al. (2002), excretion of DOC by the herbivore component is 

co = r H (A3i) 


where r H is the herbivore excretion rate at 20°C, and H 0 is the half-saturation constant for 
excretion (Table 1). H 0 is adjusted from units of pM carbon in Aumont et al. (2002) to mg m~ 3 
chlorophyll to conform to the NOBM units for herbivores. Bacterial degradation of DOC is 
represented by 


cp = A, 


DOC 


DOC 


DOC R, 


(A3 2) 


where A D oc is the DOC remineralization rate and Ki and K 2 are half-saturation constants for 
remineralization (Aumont et al., 2002; Table 1). Aumont et al. (2002) used phosphate, so here 
we substitute nitrate, since phosphate is not available in NOBM. Again parameters are allowed 
to vary as a function of temperature. In addition, the value for Ki was increased by a factor of 
10 to convert to nitrate rather than phosphate. According to Conkright et al. (1994), nitrate 
contours generally follow phosphate but nitrate concentrations are approximately 10 times 
higher. 

Dissolved inorganic carbon has a single sink, uptake by phytoplankton during photosynthesis, 
and sources deriving from respiration by phytoplankton Q in the process of growth, herbivores © 
at all times, and bacteria <p in the process of degrading DOC. There is also an interaction with 
the atmosphere (AOC 02 ) which can be a source or a sink depending upon the difference in partial 
pressures of CO 2 in the ocean and atmosphere (ApC 02 ), and the ability for gas to transfer across 
the ocean surface interface. These complex processes are follow procedures described by the 
Ocean Carbon Model Intercomparison Project (Dutay et al., 2002; Doney et al., 2004; 
http://www.ipsl.jussieu.fr/OCMIP/). 
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Table 1. Notation, parameters, and variables for NOBM. Values are provided for the parameters and ranges are 
provided for the variables. When a parameter varies according to temperature, the value at a specified temperature 
is shown and identified. Nutrient/chlorophyll ratios are variable because of photadaptation-dependence, and only 
the range is shown, corresponding to low-, and high-light adaptation, and therefore also corresponding to 
C: chlorophyll ratios of 20 to 80 g g" 1 . All time units are converted to s" 1 before use in A(l) to A(1 1). 


Symbol 

Parameter/V ariable 


Value 

Units 

General 

K 

Diffusivity 


Variable 

m 2 s" 1 

V 

3D Gradient operator 


none 

none 

V 

3D Vector velocity 


Variable 

m s' 1 

L 

Layer thickness 


Variable 

m 


Phytoplankton 

W s 

Sinking rate of phytoplankton at 3 1°C 




Diatoms 


0.75 

mcF 


Chlorophytes 


0.25 

m d' 1 


Cyanobacteria 


0.0085 

inf 1 


Coccolithophores 


0. 3-1.4 

m d' 1 

n 

Specific growth rate of phytoplankton maximum (p m ) at 20°C: 



Diatoms 


2.00 

d' 1 


Chlorophytes 


1.68 

d-‘ 


Cyanobacteria 


1.33 

d' 1 


Coccolithophores 


1.50 

d-‘ 

K 

Senescence rate 


0.05 

d-‘ 

ks 

Half-saturation for growth as i 

[unction of quota 

0.5I k 

pmol quanta m' s' 

Ex 

Total quanta (direct+diffuse) 


variable 

pmol quanta m' 2 s' 1 

R 

Temperature dependence for growth 

0.25-9.4 

none 

G 

Temperature-dependence for cyanobacteria growth 

0. 5-1.0 

none 

Ik 

Light Saturation 





Light level: Low (50) 

Medium (150) 

High (200) 



diatoms 90.0 

93.0 

184.0 

pmol quanta m' 2 s' 1 


chlorophytes 96.9 

87.0 

143.7 

pmol quanta m' s' 


cyanobacteria 65.1 

66.0 

47.0 

pmol quanta m' 2 s' 1 


coccolithophore 56.1 

71.2 

165.4 

pmol quanta m' s' 

Nutrients (N) 

bN,S,F 

Nutrient:chlorophyll ratio 





Nitrogen 


0.63 

pM (pg l' 1 )' 1 


Silica 


1.26 

pM (pg l' 1 )' 1 


iron 


0.03 

pm (pgr 1 )' 1 

8 

Nutrient excretion 





Nitrate 


0.0 

d-' 


Ammonium 


0.10 

d' 1 


Silica 


0.0 

d' 1 


Iron 


0.10 

d-' 
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kN 

Half-saturation constant nitrogen/carbon 




Diatoms 

1.0 

pM 


Chlorophytes 

0.67 

pM 


Cyanobacteria 

0.45 

pM 


Coccolithophores 

0.50 

pM 

K s 

Half-saturation constant silica 




Diatoms 

0.2 

pM 

K f 

Half-saturation constant iron 




Diatoms 

0.12 

nM 


Chlorophytes 

0.09 

nM 


Cyanobacteria 

0.08 

nM 


coccolithophores 

0.08 

nM 

0 

Iron scavenging rate 




Low iron (<0.06nM) 

2.74x1 O' 5 

d“‘ 


High iron (> Q,Q6n.M) 

1 .37x1 O' 3 

d" 1 

AFe 

Atmospheric deposition of iron 

0.03-967.0 

nmol m" 2 d" 1 

C:N 

Carbon:nitrogen ratio 

79.5 

pgl’ 1 (pM)’ 1 

C:S 

Carbon: silica ratio 

79.5 

dgl- 1 (dM)- 1 

C:Fe 

Carbon:iron ratio 

1800 

pg l' 1 (nM ) _1 

Herbivores (H) 

y 

Grazing rate maximum (y m ) at 20°C 

1.2 

d 11 

A 

Ivlev constant 

1.0 

(og rV 

r|i,r|2 

Heterotrophic loss rates 

0.1, 0.5 

d-‘ 

Rh 

Temperature-dependence for grazing 

0.75-2.7 

none 

Detritus (D) 

w d 

Sinking rate of detritus at 3 1°C 




Carbon/nitrogen detritus 

40.0 

m d' 1 


Silica detritus 

50.0 

m d" 1 


Iron detritus 

5.0 

m d' 1 

ClC,S,F 

Remineralization rate at 20°C 




Carbon/nitrate 

0.01 

d' 1 


Silica 

0.05 

d“‘ 


iron 

0.50 

d 1 

O 

Carbon:chlorophyll ratio 

Variable 

gg' 1 

Carbon (DOC, DIC) 

8 

Phytoplankton DOC exudation fraction 

0.05 

none 

IH 

Excretion rate of DOC by herbivores at 20°C 

0.05 

d- 1 

Ho 

Half-saturation constant for herbivore excretion of 
DOC 

0.14 

mg m' 3 

A,d 

Detrital breakdown rate at 20°C 

0.05 

d 1 

^DOC 

Remineralization rate of DOC to nitrate 

0.005 

d' 1 

CO 

Herbivore excretion function for DOC 

variable 

d“‘ 

n 

Phytoplankton respiration fraction 

0.05 

none 

0 

Herbivore respiration 

0.05 

d- 1 
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0 

Bacterial degradation of DOC to DIC 

variable 

d 11 

Ki 

First half-saturation constant for remineralization to 
nitrate 

3.0 

pM 

k 2 

Second half-saturation constant for 
remineralization to nitrate 

15.0 

pM 

A0cQ2 

Atmospheric-oceanic CO 2 Equilibration 

Variable 

patm 
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